mmm* 


mm 


p.  p.  unti 

.  V  u-JW ■>  * 

P  i.  SAX 
A  P.  BLAPOFOPD 
SEISMIC  DATA  LABOPATOPY 


POPCE  TECHNICAL  APPLICATIONS  CEPTEP 
Wff  */»#?«•.  O.C. 


ADVANCED  PESEAPCP  PPOJECTSABENCY 


••produced  by 

NATIONAL  TECHNICAL 
INFORMATION  SERVICE 

Springfiald.  V».  73151 


APPROVED  POP  PUBLIC  RELEASE;  DISTRIBUTION  UNLIMITED 


DISCLAIMER  NOTICE 


THIS  DOCUMENT  IS  THE  BEST 
QUALITY  AVAILABLE. 

COPY  FURNISHED  CONTAINED 
A  SIGNIFICANT  NUMBER  OF 
PAGES  WHICH  DO  NOT 
REPRODUCE  LEGIBLY. 


Mil 


GENERALIZED  LINEAR  FILTERING  OF 
SEISMIC  ARRAY  DATA 


SEISMIC  DATA  LABORATORY 

AFTAC  Project  No. : 

Project  Title: 

ARPA  Order  No . : 

ARPA  Program  Code  No. : 

Kame  of  Contractor: 


REPORT  NO.  269 

VELA  T/2706 

Seismic  Data  Laboratory 
1714 

2  F  - 1 0 

TELEDYNE  GEOTECII 


Contract  No.:  F336S7-72-C-0009 

Date  of  Contract:  01  July  1971 

Amount  of  Contract:  $  1,290,000 

Contract  Expiration  Date:  30  June  1972 

Project  Manager:  Royal  A.  Ilartenberger 

(703)  836-7647 

P.  0.  Box  334,  Alexandria,  Virginia 


APPROVED  FOR  PUBLIC  RELEASE;  DISTRIBUTION  UNLIMITED 


I  limr  I'hi 


GhNIiRALI  ZHIi  L I MliAIi  l:l  I.TliRIN'C  OF 
SI:  I  SMI  C  A  RIMY  DATA 


•  wHwncznszrtsnssr 


l'1,ltZ»  I),R*i  Sax»  K.L.J  and  Bland  ford,  R.R. 


I  WWHTBtTI 

October  I97i 

i#  CONTRACT  04  OIIANT 

F33GS7-70-C-0941 

VI;  LA  T/'  07  06 

«. 

ARPA  Order  No.  624 
-■  ARPA  Program  Code  No,  9 mo 


■■ 


I « a  ■»  1-0  r>n  m  n  t  *ju 


[**  _ _ _ 


APPROVED  for  PUBLIC  RELEASE.  OISTRIRUTIOR  UR  LIMIT  ED 


I II  A || T MAC T 


^‘v“"",T4!s"L)"m'i)().ii;(:ts  acp 


*  Tltrr  irncrolttp 
Ihy  il«pl£  hPijnr  i , 


•«  KEY  WO  AOS 

Array  Processing 
l-ogari  thin ic  Deconvolution 


(lone  rn  1  i  zed  Linear  Filtering 
lloinomorph  ic  Filtering 


WNCLASSirilil) 

Security  Cleanfication 


AifiTMCf 


Ifc»c  ifitory  of  f**#r*l is»4  fHs*rl«*  f0j9*ofc*feo* 

1966,  1969)  U  aryl i*4  *•  \h*  yrol*l*w  of  «*r*|iii  ih* 

xramf «r  fuat'ioft*  *  «*»***«  **«*«*  »**  ******* 

r*cor4io|  »talio^  «•  a  eo*tirMrMai»*f£*4  array# 

To  pr«ri4e  ftvoihoiic  4aia  fer  a  «*»«  ®f  »**  yroet§*or. 
a  *e«j*ic  *«ur«*  *!*»»*  »*  J**»**4  through  a to*  differed 
ra«4om  ye  normal  ion*  of  a  **locity-4*rifc  ai»«c*»r*  yicUins 
wlnr  different  *y«»6*uc  U&  owcUar  r#plo*fe»* 

fLOKC  $M9T<  »«4  tttUMMfi)  a*s4  *»  tftJrcawof f  i*l»»4t.  **‘*o* 

(2:  Nevcfttwr  19**1  6f«  *1*®  procr***®  «*»»*  *****®**  *r°* 

tbe  fc.oo«-»-»i!»*r-SrMflsU»Wc*594»fflrRi  ©#t.  *«**  of  ****  r*,**’r“ 

fetralieft*  art  f*oo«4  Ifc*  *>©56*«»«  **i**®$»*«  hf 

t«n«ftali i*4  I linear  fittrti**.  1fc«  **F»®*»®*  **4  the  **rth- 
^al»  *e»*oo#ra**  *f*  *iWf»ifl*4  kf  tM* 

•ffco  e«ftti<|uaife  to  k**'*  iow<t*r  I***  *<’p**  r*tMtft 

i**  jfirff  gcwrraiUrt*  lift*#*  filter i**  *»4  ,fcr  mw<l*4f  rt'rM 
,,m*r  «»r*  *  i  ni  Ur  feeUtr  ifcc  *r«l«'«l  «f  l>»*  ^ilH 


TABLE  OF  CONTENTS 


Page  No. 

ABSTRACT 

INTRODUCTION  1 

Seismic  signal  models  1 

Tine  probability  density  of  the  random 

transfer  functions  under  the  hypothesis 

of  zero  additive  noise  2 

The  optimum  processor  5 

RESULTS  7 

Simulated  receiver  site  effects  7 

Observed  data  9 


REFERENCES 


12 


LIST  UF  FIGURES 


Figure  Title 

Transformation  of  a  normal  probability  density  into  an  almost 
uniform  probability  density. 

The  best  fit  to  a  normal  probability  density  of  Klappenberger 's 
(19o7)  data.  (After  Klappenberger,  1967). 

Klappenberger' s  best  fit  of  LASA  peak-to-peak  amplitude  data 
to  a  cumulative  normal  distribution  function. 


Figure  No. 
1 

2 

3 


A  schematic  representation  of  the  normal  incident  P-wave  layering 
problem. 

Analog  computer  simulation  of  the  transmission  problem.  Note 
the  similarity  of  this  problem  to  the  parametric  amplifier 
problem.  Also  note  that  S-waves  could  be  handled  by  designing 
another  ladder  network  with  different  delays  and  reflection 
coefficients,  and  interconnecting  the  nodes  of  both  networks 
with  the  proper  conversion  factors. 

A  two-stage  generalized-linear  filter  for  the  estimation  of  a 
signal  convolved  with  equal  variance,  log-normally  distributed, 
uncorrelated  transfer  functions. 


A  modification  of  st^ge  2  of  the  GLF  in  Figure  6  to  delineate 
long  or  short  time  components.  Tapering  the  complex  cepstrum  is 
quite  often  done  in  speach  processing  (Noll,  1964;  Oppenheim,  1968). 

Synthetic  seismogram  1  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs.  depth  (c)  Velocity  vs.  depth. 

Synthetic  seismogram  2  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs.  depth  (c)  Velocity  vs.  depth. 

Synthetic  seismogram  3  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs.  depth  (c)  Velocity  vs.  depth. 

Synthetic  seismogram  4  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs .  depth  (c)  Velocity  vs.  depth. 

Synthetic  seismogram  5  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs .  depth  (c)  Velocity  vs.  depth. 

Synthetic  seismogram  6  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs.  depth  (c)  Velocity  vs.  depth. 

Synthetic  seismogram  7  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs.  depth  (c)  Velocity  vs.  depth. 


LIST  OF  FI  CUKES  (Cont'tl.) 


Figurt  Title 

Synthetic  seismogram  8  (a)  Reflection  coefficients  vs,  depth 

(b)  Seismogram  vs.  depth  (c)  Velocity  vs.  depth. 

Synthetic  seismogram  9  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs.  depth  (c)  Velocity  vs.  depth. 

The  nine  synthetic  seismograms  and  the  source  signal  vs.  time. 

Synthetic  seismogram  1  (a)  Seismogram  ys.  time  (b)  Log-amplitude 

vs.  frequency  (c)  Unwound  phase  vs.  frequency 

Synthetic  seismogram  2  (a)  Seismogram  vs.  time  (b)  Log-amplitude 

vs.  frequency  (c)  Unwound  phase  vs.  frequency. 

Synthetic  seismogram  3  (a)  Seismogram  vs.  time  (b)  Log-amplitude 

vs.  frequency  (c)  Unwound  phase  vs.  frequency. 

Synthetic  seismogram  4  (a)  Seismogram  vs .  time  (b)  Log-amplitude 

vs.  frequency  (c)  Unwound  phase  vs.  fi^qucncy. 

Synthetic  seismogram  5  (a)  Seismogram  vs.  time  (b)  Log-amplitude 

vs.  frequency  (c)  Unwound  phase  vs.  frequency. 

Synthetic  seismogram  6  (a)  Seismogram  vs .  time  (b)  Log-amplitude 

vs.  frequency  (c)  Unwound  phase  vs.  frequency. 

Synthetic  seismogram  7  (a)  Seismogram  vs.  time  (b)  Log-amplitude 

vs.  frequency  (c)  Unwound  phase  vs.  frequency. 

Synthetic  seismogram  8  (a)  Seismogram  vs.  time  (b)  Log-amplitude 

vs.  frequency  (c)  Unwound  phase  vs.  frequency. 

Synthetic  seismogram  9  (a)  Seismogram  vs.  time  (b)  Log-amplitude 

vs.  frequency  (c)  Unwound  phase  vs.  frequency. 


Figure  No. 


Results  of  processing  the  nine  synthetic  seismograms 

(a)  Weighted  beamforming 

1.  Seismogram  vs.  time 

2.  Log-amplitude  vs.  frequency 

3.  Unwound  phase  vs.  frequency 

(b)  Generalized  linear  filter  1 

1.  Seismogram  vs.  time 

2.  Log-amplitude  vs.  frequency 

3.  Unwound  phase  vs.  frequency 

4.  Complex  cepstrum  vs.  time 

(c)  Generalized  linear  filter  2 

1.  Veismograin  vs.  time 

2.  Log-amplitude  vs.  frequency 

3.  Unwound  phase  vs.  frequency 

4.  Tapered  complex  cepstrum 


LIST  OF  Flames  (Cool'd.) 


Fi£uro  So. 

:s 


Figure  Title 

Location  of  events  and  stations. 

LONGSIIOT,  Station  KX-UT  (a)  Seisaogran  vs.  tine,  (b)  Log- amplitude 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

LONGSIIDT,  Station  LC-XM  (a)  ScisBograxs  vs.  tine,  (b)  Los-asplitudo 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

LOXGSIIOT,  Station  SV3QB  (a)  Seismogram  vs.  time,  (b)  Log-amplitude 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

LOXGSIIOT,  Station  IIL2ID  (a)  Seismogram  vs.  time,  (b)  Log-amplitude 
vs.  Itirequcncy,  (c)  Unwound  phase  vs.  frequency. 

LOXGSIIOT,  Station  UG-SU  (a)  Seismogram  vs.  time,  (b)  Log-amplitude 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

LOXGSIIOT,  Station  UK-ON  (a)  Seismogram  vs.  time,  (b)  Log-amplitude 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

LOXGSIIOT,  Station  SJ-TX  (a)  Seismogram  vs.  time,  (b)  Log-amplitude 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

LOXGSIIOT,  Station  TF-CL  (a)  Seismogram  vs.  time,  (b)  Lop-ampl ituJe 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

LOXGSIIOT,  Station  KC-MO  (a)  Seismogram  vs.  time,  (b)  Log-amplitude 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

LONGSIIOT,  Station  YR-CL  (a)  Seismogram  vs.  time,  (b)  Log-amplitude 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

LONGSIIOT  processing  results 

(a)  Beamforming 

1.  Seismogram  vs.  time 

2.  Log-amplitude  vs.  frequency 

3.  Unwound  phase  vs.  frequency 

(b)  Generalized  linear  filter  1 

1.  Seismogram  vs.  time 

2.  Log-amplitude  vs.  frequency 

3.  Unwound  phase  vs.  frequency 

(c)  Generalized  linear  filter  2 

1.  Seismogram  vs .  time 

2.  Log-amplitude  vs.  frequency 

3.  Unwound  phase  vs.  frequency. 

M1LR0W,  Station  WQ-IL  (a)  Seismogram  vs.  time,  (b)  Log-amplitude 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 


20 

30 

31 

32 

33 

34 

35 

36 

37 

38 


39 

40 


US!  Of  9WM$  fCtou a^) 


I  tigwc  I  Hie 

MlUdf,  SliiiM  It- 1©  (a)  S*i&Mtr»¥»,  timr.  (%»| 
v*.  fftfiKstj1,  (e)  ito*34»>l  fhM-G  v* ,  frifmiy 

HILMW.  Station  $J*1*  (a)  5<fl*m»jgf.wt  «*,  umr,  (%)  tQf»a*e lliwlt 
v*.  frt^KRqr,  (c)  Ow^nl  ^u«  w,  frcf^y 

uiusotf,  Su: ica  WiJii  (a)  Seiswafrain  v*.  tin#,  |*>J 
v*.  frequency,  (cj  phase  *■»,  fre-fs *t*«y 

MiUttflJf,  Station  IC-H'i  (a)  Seismogram  time,  |bj  log  amytimtie 
v*.  frequency,  (tl)  Unwound  phaae  « ,  IrtpHKy. 

HIIJUjM,  Station  KN-Uf  (a)  Seismogram  v*.  tine,  |b)  t«g>a*f Siltfefc 
vs.  frequency,  (c)  Unwound  phase  v*»  frequency. 

’IILHWt*,  Station  AS-PA  (a)  Seismogratt  v*.  time,  (b)  l«fa«filiia4r 
vs.  frequency,  (c)  Unwound  phase  v*.  frequency. 

Ml I.RUft,  Station  LU2AL  (a)  Seismogram  v*.  time,  (b»  Cng-amptluale 
Vi.  frequency,  (c)  Unbound  phase  vs.  frequency. 

MILKOK,  Station  PJ-PA  (a)  Seisrogran  vs.  Use,  (b)  log -amplitude 
vs.  frequency,  (c)  Unwound  phase  vs.  frequency. 

Processing  results  for  MILKW 

(a)  Ucanforming 

1.  Seismogram  vs.  tine 

2.  Log-anpl itude  vs.  frequency 

3.  Unwound  phase  vs.  frequency 

(b)  Generalised  linear  filter  ! 

1.  Seismogram  vs.  tine 

2.  Log-amplitude  vs.  frequency 

3.  Unwound  phase  vs.  frequency 

4.  Complex  cepstrum  vs.  time 

(c)  Generalised  linear  filter  2 

1.  Seismogram  vs.  time 

2.  Log -amplitude  vs.  frequency 

3.  Unwound  phtnc  vs.  frequency 

4.  Tapered  complex  cepstrum. 

22  Nov.  1965  Andreanoff  Islands  event,  Station  liN-MO 
(a)  Seismogram  vs.  time,  (b)  Log-amplitude  vs.  frequency, 

(c)  Unwound  phase  vs.  frequency. 

22  Nov.  1965  Andreanoff  Islands  event,  Station  IIV-MA 
(a)  Seismogram  vs.  time,  (b)  Log-amplitude  vs.  frequency, 

(c)  Unwound  phase  vs.  frequency. 


U£V  m  W9M&  |Cmi"4.| 


ft  i*8lf 

22  <*«*•  W«fcV  |*|M«t:fc  <miH. 

6*1  Hn  iiflut,  tuf^'SUiilt  i*»  ft tmmm  t, 

If  9  MwwNfti  fta*#  m»  fit 

22  %m  >  I  Mil  ,»A#rw»##f  lil.iithl*  $9*5  lift  4ft  WS 

1*9  Sf  l****#f*i»  vv*.  iH*t,  !«<•  MflHrtti*'  rt  ftwmew'ft 
|c|  m*®m$  jru**  »•*  ftvyvtmif 

22  %m*  t'Hr  r«,  -i.ium  *»*■■%* 

1*9  **#  ihhp,  |9»|  *•#  Uwmmi** 

tel  f#1**#  v»*.»  fl'StfMJiWf"  ■ 

22  .V.'*  m*  A%fel'5f!W«f9  i*f*mt*  «w*!| ,  $9*9 !•*  Ht4-3> 

it  Hi#  !*<rll*>mff9;M*  4*  J  IMr ,  $|>f  fiteif  vi'* 

$4$  WfrWBflj  fftkMt-  4*» 

22  ftsft-  rffc*  <teJt«M*ff  Mil  i«  Miftt 

1*9  «**  im«,  I9"9  tat  amflimnif  «■* 

fCf  ifUMmiM  119***  4'* 

22  ta*  91*4  9*l4rtt#*  **«*»?.  *l*9»*»«  Wft 

1*9  »'fi  **«**■?>**'  vr*  I  Hr,  l&i  l«j|-  .Wif9ll'«#e  *'*  0«*9N«*fy, 

Ik 9  swwtftHl  |6<«e  #*»  ti&y MMf* 

9"#«*r  *  *  #*•*••*$«*  Hi 

1*9 

I,  StlHWp^fr  «*  iiifw*' 

*.  lifiiMk1  vt'*a  tH'tf'H'Mii 

i.  fi***  v*  ft’t'vwmf 

1*9  ttHKffttilaJ  I9»**f  9 

la  **  l*!»? 

«•  llkf'llg  *|tflf8*lte  *'* 

I»  Utow&auS  |kfc*<*'  v* 

>i.  twmfij**  ;  tw  r*~  iM* 

k|  ft***?  filler  ? 

la  r*  I  Hi* 

c«  l*f  «'*  f  f(8sfii¥M* 

i*  hw»*tft3  ffi.**-*  »•*« 
i  I*f*i*tJ  t*my9f  »■.  4*^*9  fm». 


iQSWff  ta. 


it  * 

&r 


S3 

$1 


*a 


a-*/ 

Jv  T 


i* 


its*  or  v«Mf$ 


lisi* 

oh#!  it,  tM» 


li*VI<e-  Vt>. 

I 

n 


INTRODUCTION 


Seismic  signal  models 

One  difficulty  encountered  in  beaming  a  short-period 
worldwide  or  cent inenta 1 -s i zed  seismic  array  is  that  the 
signal  coherency  falls  off  markedly  with  distance* 
(Hartenberger,  1967;  K1 appenberge r,  1967;  Capon  ct  a  1 , 

1967). 

A  possible  model  of  the  seismogram  recorded  at  the 
jth  instrument  an  array  under  excitation  of  a  signal  plus 
uncorrclatcd  additive  noise  is 

Yj(t)  a  s(t J»hj (t)  ♦  n j (t) 

where  the  *  indicates  convolution. 

The  Fourier  transform  vf  (1)  yields 

Yj(f)  =  S(f)llj(f)  ♦  Nj(0 

The  Hj ( f )  may  be  thought  of  as  a  source-receiver  travel-path 
transfer  function  for  the  jth  receiver. 


•Amplitude  anomalies  may  run  as  high  as  cight-to-onc  for 
LASA  (Klappcnbcrger,  1967;  Mack,  1969),  whereas  time 
anomalies  may  run  as  high  as  1.5  sec  at  LASA  (Chiburis, 
1966)  . 
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The  probability  density  of  the  random  transfer  functions 
under  the  J]_^[£iO|tji^s^s_ojf^ze_ro_ja^djdi^tj. ve  noise _ 

If  we  make  the  assumption  of  large  signal-to-additive 
noise  ratio,  equation  (2)  becomes 

Yj(f)  3  Nj(f)  S(f) 


If  we  further  assume  that  IT  (f)  may  be  written  as  a 
product  of  a  large  number  of  small  disturbances,  then 


Yj  (f ) 


M 

n 

m=l 


Amj(f)  I  exP  i*mj(f)]  S(f) 


(4) 


where  m  is  the  index  of  each  small  disturbance  and  M  is  the 
total  number  of  disturbances.  Taking  the  logarithm  of  both 
sides  of  equation  (4),  gives 

M 

log  Yj  (f)  =  log  |  S  (f )  |  +  £  logl^ff)! 

m=l  J 

*  1  |*s<«  *  ,£  [♦„,(«  ±  2  7,Mmj(f)j| 

where  Mm..(f)  is  an  integer  which  reflects  the  uncertainty  in 
the  correct  Riemann  sheet  at  each  frequency.  We  see  that  the 
logarithm  of  the  complex  Fourier  spectrum  has  converted  our 
multiplicative  variables  into  a  new  set  of  additive  variables. 

If  Amj(^f-)  and  the  ( f )  can  be  described  by  a  probability 
density  function,  then  we  may  invoke  the  weak  law  of  large 
numbers  and  state  that  our  new  random  variable  log  Y.(f) 
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is  asymptotically  normally  distributed  over  the  ensemble 
with  the  expected  value  of  the  real  part  equal  to 

log  I S ',f)  I  ♦  log|u(f)  I 

where  p  ( f )  is  the  geometric  mean  of  the  random  transfer 
functions.  The  expected  value  of  the  imaginary  part  of  log 
Y j (f)  is  equal  to 

‘f’s  (f)  ♦ 


where  <J>m(f)  is  the  phase  of  the  iTean  of  the  random  transfer 
functions.  The  maximum  likelihood  estimate  of  log | p ( f ) |  is 
given  by 


N 


log  |  p  (f )  I  =  L 


M 

£ 


*  jti  A  l0E'Vf)l 


(6a) 


The  maximum  likelihood 


estimate 


of  4>u(f)  is 


given  by 


1  i-m  '  i  | 


M 

£ 

m=l 


^mj  ^ 


(6b) 


The  probability  density  function  which  describes  the 
complex  random  variable  y  ^  is  a  log-normal  density  function 
(Aitchinson  and  Brown,  1963).  Note  that  if  the  signal 
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spectrum  is  white,  then  the  estimates  at  different  values 
of  frequency  are  normally  distributed  and  we  have  a 
stationary  bivariate  log-normal  variable.  If  the  signal 
spectrum  is  non-white  we  have  a  non-frequency-stationary 
bivariate  log-normal  variable. 

The  theory  of  logarithmic  deconvolution  has  been 

developed  by  Oppenheim  (1966)  and  Schafer  (1967,  1969)  as 

a  specific  case  of  Oppenheim's  generalized  linear  filtering 

theory.  For  random  transfer  functions,  the  real  difficulty 

with  this  method  is  the  choice  of  the  proper  Riemann  number 

Nj(f)  where  N^(f)  =  £  Mmj  (f)  for  each  recording  site.  If  the 

m=  1 

phase  changes  rapidly  with  frequency,  i.e.,  ( f )  has  a  large 

variance,  we  may  find  ourselves  faced  with  an  almost  uniform 
distribution  (over  the  ensemble  of  stations)  for  our  estimate 
of  the  ensemble  average  of  0(f)  due  to  the  discontinuous 
Riemann  transformation.  This  problem  is  illustrated  in  Figure 
1  for  a  zero  mean  distribution.  Our  solution  is  to  set  the 
phase  equal  to  zero  zero  at  zero  frequency  for  each  recording 
site  and  then  unwind  the  phase.  Schafer  (1969)  has  discussed 
techniques  for  unwinding  the  phase. 

There  is  experimental  evidence  that  the  A.(f)  are  log- 
normally  distributed,  Klappenberger  (1967)  has  shown  that 
earthquakes  have  log-normally  distributed  peak-to-peak  ampli¬ 
tudes  at  LASA  and  Freedman  (1967)  has  shown  the  same  distri¬ 
bution  for  earthquakes  at  worldwide  recording  sites.  Figure  2 
is  the  histogram  of  K1 appenbe r ger ' s  data.  Figure  3  is 
Klappenberger' s  best  straight  line  fit  to  his  data  plotted 
on  log-normal  probability  paper.  If  the  data  are  log-normal, 
the  cumulative  frequency  distribution  should  approximate  a 
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straight  line.  The  hypothesis  that  the  data  are  log-normal 
cannot  be  rejected  at  the  95%  confidence  level. 

Further  indication  of  the  log-normality  of  the  station 
spectra  may  be  obtained  by  passing  a  synthetic  signal 
through  ensembles  of  plane  parallel-layered  media  with 
random  reflection  coefficients. 

Figure  4  illustrates  the  problem  of  a  seismic  wave 
impinging  upon  a  set  of  plain,  parallel  layered  media.  A 
subroutine  given  by  Claerbout.  (1968)  was  programmed  to 
yield  the  transmission  seismogram  of  an  arbitrary  wavelet 
impinging  upon  the  bottom  layer.  Figure  5  is  the  analog 
computer  simulation  of  the  transmission  and  reflection 
problem.  Electrical  engineers  might  recognize  the  circuit 
as  a  "ladder  network"  or  as  a  parametric  amplifier. 

The  optimum  processor 

Under  the  assumption  of  white,  uncorrc 1 ated ,  zero-mean, 
equal  variance,  log-transfer  functions  and  zero  additive  noise 
it  can  be  demonstrated  from  maximization  of  the  1 og- 1  ike  1 ihood 
function  that  the  likelihood  estimator  for  the  logarithm  of 
the  signal  in  the  frequency  domain  is  just  the  arithmetic 
average  of  the  individual  station  log-spectra. 

^  N 

log  S(f)  »  i  £  log  Y.(f)  (7) 

N  j  =  1  J 

or  S(f)  =  exp  £log  S(f)J  (8) 
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O) 


_  "  N  "Wn 

so  that  S(f)  =  ||  Y ;  (f) 

3=1 


Equation  (9)  is  just  the  geometric  mean  of  our  individual 
station  amplitude  spectri.  Figure  6  is  a  block  diagram  of 
the  processor  for  this  case.  Figure  7  is  a  modification  of 
the  processor  to  filter  out  undesirable  time  components. 
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beat?  summing.  Trace 


/  N  /N  v 

al  is  the  weighted  beam  sum  l/L  Wj  Xj  (t)/2Z  W;  }  I  where 
*  l/o^  and  where  is  the  sample  standard  deviation. 

Trace  a2  is  the  log  amplitude  of  the  weighted  beam  sum 
versus  frequency.  Trace  a  3  is  the  phase  of  the  weighted 
beam  sum  versus  frequency.  Trace  b!  is  the  geometric  mean; 
trace  b2  is  the  averaged  log-amplitudes  versus  frequency; 
and  trace  b3  is  the  average  phases.  Trace  b4  is  the  complex 
ccpstrum  (Schafer,  1369) .  Trace  b4  appears  to  be  anti¬ 
symmetric  about  the  midpoint  because  the  phase  (trace  b3) 
is  much  larger  than  the  log  amplitude  (trace  b4).  The  unwound 
phase  is  the  imaginary  part  of  the  log-spectra  and  is  anti¬ 
symmetric  about  its  midpoint  (10  Hz),  and  the  log  amplitude 
Is  the  real  part  of  the  log  spectra  and  is  symmetric  about 
its  midpoint  (10  Hz).  Since  the  numbers  for  the  unwound 
phase  are  much  larger  than  the  numbers  for  the  log-amplitude, 
the  complex  ccpstrum  appears  to  be  antisymmetric. 


The  traces  under  c  have  been  processed  by  tapering  in 
the  pseudo-time  domain.  Trace  c2  is  the  smoothed  averaged 
log-amplitudes,  and  trace  c3  is  the  average  unwound  phase. 
Trace  c4  is  the  first  half  of  the  tapered  complex  cepstrum. 
It  would  appear  that  the  generalized  filters  have  done  a 
better  job  of  estimating  the  signal  than  has  the  weighted 
beamed  sum.  The  log  amplitudes  of  the  geometric  beams  and 
the  phases  of  the  geometric  beams  also  show  less  variance 
than  do  the  log-amplitude  and  phase  of  the  weighted  beam 
sum.  Notice,  hoittvef,  the  small  precursor  introduced  by  the 
generalized  fillters. 
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Observed  data 


Three  events  were  processed  on  a  continental  sized  array 
of  LRSM  (Long-Range-Seismic-Measurement)  stations  distributed 
throughout  the  continental  United  States  and  Canada,  Two 
events  (LONG  SHOT  and  MILROW)  were  underground  nuclear  explo¬ 
sions  which  were  detonated  on  Amchitka  Island  in  the  Aleutian 
Island  Arc,  The  other  event  was  an  Andreanoff  earthquake 
(22  November  1965)  in  the  same  island  arc.  Figure  28  is  a 
map  of  the  three  events  and  the  recording  stations.  Table  I 
gives  the  seismic  parameters  of  interest  for  the  three  events. 
Table  II  gives  the  locations  of  the  recording  stations  used 
in  this  study.  The  processing  results  for  LONG  SHOT  are 
given  in  Figure  29  through  Figure  37.  The  formats  are  the 
same  as  for  Figures  18  through  27.  Figures  29  through  37  are 
for  stations  KN-UT,  LC-NM,  SV3QB ,  HL2ID,  RK-ON,  SJ-TX,  TF-CL, 
KC-MO,  and  YR-CL,  respectively.  Notice  the  large  variations 
in  signal  waveform  from  station-to-station  in  all  traces. 

In  Figure  39,  note  the  small  precursor  introduced  by  the 
processor  since  we  did  not  demand  minimum  phase  as  one  of 
the  attributes  of  our  processor.  The  main  signal  of  our 
processor  (trace  39bl)  seems  to  "ring"  more  than  the  signal 
form  the  weighted  beam  (trace  39al),  yet  the  coda  dies  out 
more  rapidly  for  the  first  generalized  linear  filter  compared 
to  the  weighted  beam.  The  variance  of  trace  b2  and  trace  b3 
is  seen  to  much  less  than  the  variance  of  trace  a2  and  a3. 

The  spectral  null  at  1.8  Hz  in  the  LONG  SHOT  log-  i  -e 

amplitude  spectra  (traces  a2,  b2,  and  c2)  has  been  inter¬ 
preted  (Cohen,  1969)  as  interference  between  the  direct  P 
wave  and  pP,  the  reflection  from  the  free  surface. 
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The  processing  results  for  the  MILROW  event  are  given 
in  Figures  40  through  43.  In  Figure  40  through  Figure  48, 
we  have  the  results  (in  the  same  format  as  for  the  LONG 
SHOT  event)  for  stations  WQ-IL,  BY-IO,  SJ-TX,  CR2NB,  LC-NM, 
KN-UT,  AS-PA,  EU2AL,  and  PJ-PA.  Figure  49  contains  the 
results  of  our  three  processors  in  the  same  format  as 
before.  Notice  that  the  first  spectral  null  is  shifted  to 
a  lower  frequency  from  the  LONG  SHOT  null  (trace  37c2).  This 
is  consistent  with  the  longer  time  delay  for  the  reflected 
phase  due  to  the  deeper  depth  of  MILROW  compared  to  LONG 
SHOT.  Notice  the  striking  similarity  between  the  processed 
LONG  SHOT  event  (trace  37cl)  and  the  processed  MILROW 
event  (trace  49cl)  despite  the  great  variablity  of  the 
individual  stations.  The  decay  times  appear  to  be  identical. 
The  difference  in  the  two  waveforms  appears  to  arise  from 
the  pP  phase,  since  the  first  half  cycle  of  motion  (neglecting 
the  precursor)  is  similar  in  shape  and  frequency.  The  results 
for  the  22  November  1965  Andreanoff  earthquake  are  contained 
in  Figure  50  through  Figure  58.  Figures  48  through  55  are 
for  stations  EN-MO,  HV-MA,  KC-MO,  KN-UT,  MN-NV,  RG-SD,  RK-ON, 
and  WN-SD. 

Comparing  the  result  for  the  Andreanoff  Islands  event 
(trace  58cl)  with  the  MILROW  event  (trace  49cl)  and  the  LONG 
SHOT  event  (trace  37cl),  we  see  that  the  earthquake  signal 
is  more  complex  and  exhibits  a  longer  decay  time  than  either 
MILROW  or  LONG  SHOT.  However,  this  fact  is  also  apparent  at 
the  individual  stations  and  on  the  weighted  beams. 

In  conclusion  we  state  that 

The  generalized  linear  filters  did  at  least  as  well  as 
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the  weighted  beam  sum  in  the  time  domain  in  estimating  our 
synthetic  signal.  No  significant  signal  distortion  in  the 
first  motion  was  present  in  our  estimate. 

The  generalized  linear  filter  traces  for  the  explo¬ 
sions  were  similar  despite  substantial  variations  in  the 
records  for  the  same  event  at  different  stations.  The 
explosion  traces  were  simpler  than  the  earthquake  trace, 
but  this  was  also  true  at  the  individual  stations. 

The  codas  are  smaller  for  the  output  from  the 
generalized  filters  than  at  the  individual  stations  or  on 
the  weighted  beam  sums.  The  log-amplitude  and  phase  spectra 
of  the  generalized  linear  filter  output  exhibited  less 
variance  from  their  mean  compared  to  either  individual 
station  spectra  or  the  weighted  beam  spectra.  This  might 
be  expected  since  the  geometric  beam  is  the  minimum  vari¬ 
ance  estimator  for  uncorrelated  equal  variance  random  trans¬ 
fer  functions. 

The  results  suggest  that  continental  beamforming  or 
indeed  even  world-wide  beamforming  is  possible  for  the 
initial  P-wave  using  the  theory  of  generalized  linear 
filtering,  and  that  complexity  measurements  and  analysis 
of  source  dynamics  might  better  be  performed  on  the 
generalized  beam. 
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Stations  Used  in  this  Study 

Station 

Location 
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Figure  3.  Klappenberger's  best  fit  of  LASA  peak-to-peak  amplitude  data 
to  a  cumulative  normal  distribution  function. 


Figure  5 


Analog  computer  simulation  of  the  transmission  problem. 

Note  the  similarity  of  this  problem  to  the  parametric 
amplifier  problem.  Also  note  that  S-waves  could  be 
handled  by  designing  another  ladder  network  with  different 
elays  and  reflection  coefficients,  and  interconnecting  the 
nodes  of  both  networks  with  the  proper  conversion  factors. 
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short  time  components.  Tapering  the  complex  cepstrum  is  quite  often 
done  in  speech  processing  (Noll,  1964;  Oppenheim,  196S) . 
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Figure  10.  Synthetic  seismogram  3  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs.  depth  (c)  Velocity  vs.  depth. 
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Figure  11.  Synthetic  seismogram  4  (a)  Reflection  coefficients  vs.  depth 

(b)  Seismogram  vs.  depth  (c)  Velocity  vs.  depth. 
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Figure  5-1.  LONGSIiOT,  Station  RK-OX  (a)  Seismogram  vs.  time,  (b)  Log -amplitude  vs.  frequency 
(c)  Unwound  phase  vs.  frequency. 
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